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_ , Equilibrium vortex configuration in conventional type II superconductors containing columnar 

_ ■ defects is theoretically investigated. Near the upper critical field a single defect causes a strong local 

deformation of the vortex lattice. This deformation has C3 or Ce point symmetry, whose character 
strongly depends on the vortex-defect interaction. If the interaction is attractive, the vortices can 
collapse onto defect, while in the case of repulsion the regions free of vortices appear near a defect. 
Increasing the applied magnetic field results in an abrupt change of the configuration of vortices 
giving rise to reentering transitions between configurations with C3 or Ce symmetry. In the case 
' of a small concentration of defects these transitions manifest themselves as jumps of magnetization 

<**] , and discontinuties of the magnetic susceptibility. 

o ■ 

O ; PACS: 74.60.Ge; 74.62.Dh 

U • 
Of 

^ ' I. INTRODUCTION 

xn \ 

Mixed state or Shubnikov phase! of type II superconductors is characterized by penetration of 
vortices into the sample! each one carrying the superconducting flux quantum <ps = vr/ic/e. A single 
vortex has a normal core with radius of order £(T) (the coherence length at temperature T) sur- 
rounded by a closed superconducting current occupying the tube with radius A(T) (the penetration 
length). Near the upper critical field H C 2 these vortices form a triangular Abrikosov lattice!. If an 
external current is applied, vortices start to move due to the Lorentz force. This motion leads to 
an energy dissipation. Different kinds of defects, such as dislocations, point defects or regions with 
different superconducting properties create some additional field acting on the vortices. As a result 
vortices are pinned and nondissipative current of finite amplitude can flow through the superconduct- 
ing sample. There are two different types of pinning!. Weak defects lead to the so-called collective 
pinning^. In this^regime the vortex lattice is slightly deformed. This deformation is well described by 
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elasticity theoryUm Strong defects lead to a single-particle pinning%. A single strong defect is able 



to pin a vortex and, at a finite defect concentration, formation of metastable states is possibles. A 
detailed theory of vortex pinning in conventional type II superconductors was formulated by Larkin 
and Ovchinnikov (see review paper!). 

The discovery of high-temperature superconductivity^ resulted in a deep understanding of the new 
rich and fundamental properties of vortex systems (see an exhaustive review papers of Blatter et. 
a/.ty and BrandtEJ and references therein). Statistical mechanics of vortices was formulated and 
new concepts appeared such as melting of the vortex lattice, vortex liquid and vortex glass. The 
^ ■ usage of heavy ion irradiation for preparation of superconducting samples with columnar defects^ 
opened new experimental possibilities to study the properties of vortex matter. Columnar defects 
serve as strong pinning centers, each of which is able to pin a single vortex as a whole. The radius 
of the columnar defect could be less than the Abrikosov lattice constant a near H C 2. Such defects 
are referred to as short-range ones. 

Strong columnar defects with radius L much larger than the coherence length may lead to the 
formation of multiquantum vortices in high temperature superconductors!!!. Such vortices were ob- 
served experimentally on submicron artificial holes in mutlilayers Pb/Ge^. Multiquantum vortices 
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can also be formed at large pinning centers with radius of order of the penetration lengthlli In 
this paper we show that columnar defects can also strongly affect the properties of conventional 
type II superconductors. In such superconductors near the upper critical field even the short-range 
columnar defects cause a strong local deformation of the vortex lattice due to its softening and as a 
result to the formation of multiquantum vortices. 

In the main part of the paper we consider a superconductor containing a single short-range colum- 
nar defect. When the magnetic field approaches to the strength of the defect effectively increases 
resulting in strong lattice deformation in its vicinity. Initially the Abrikosov lattice is triangular. 
Therefore the local deformation belongs to one of the two possible symmetry types - C$ or C3. 
In the case of an attractive defect, the vortices can collapse onto this defect with increasing of a 
magnetic field. As a result, reentering transitions between two local symmetries are possible. For 
example, at some external field the local vortex configuration Cq with a single vortex pinned by a 
defect is preferred over the C3 configuration with the defect placed at the center of a triangle. But 
with increasing of a magnetic field the closest three vortices in the C3 configuration could collapse 
onto the defect and a C^—C^ transition occurs. Now the mostly preferred local configuration is of C3 
type, with a three-quantum vortex on the defect. Further increasing of a magnetic field results in a 
C*6 type configuration with a seven-quantum vortex at the defect and so on. In the case of a small 
concentration of defects (which was realized in an experimenilll, where the radius of the defect is 
equal to 2.5 nm and the average distance between defects is 4,600 nm). These transitions manifest 
themselves as jumps of magnetization and discontinuties of the magnetic susceptibility curve. 

The present paper is organized as follows. In the second section we formulate the problem. Fur- 
ther, in the third section we study the relatively simple case of small deformation of the vortex 
lattice. It is realized for weak defects or for values of the magnetic field which are not very close to 
H c2 . The results of this study enable us 1) to scale the defect parameters with the magnetic field, 
and 2) to predict the occurrence of symmetry change when the applied magnetic field increases. The 
central section IV containes the results of the numerical solution of the pertinent equations. Here we 
present a universal phase diagram of the superconductor near the upper critical field, and analyze 
the vortex lattice deformation as a function of the magnetic field for various defect parameters. 
The next section V is devoted to the case of small concentration of the defects. Here we estimate 
the high order concentration corrections with respect to defect and study the magnetization and 
the magnetic susceptibility behavior near the upper critical field. Section VI summaries the main 
results. In the Appendix, the Abrikosov lattice expansion in terms of the first Landau level wave 
functions is obtained for an arbitrary position of the lattice with respect to the origin. 



Consider a superconductor containing columnar defects and subject to an external magnetic field 
H = Hz. Both the defect column axis and the magnetic field are assumed to be directed along the 
z-axis. The unit volume thermodynamic potential of such a superconductor at a fixed temperature 
T close to the critical temperature T c can be written as 



where r is the two dimensional (2D) position vector and the Ginzburg-Landau density / of the 



II. FORMULATION OF THE PROBLEM 
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/ (r; [a, 7, ¥, A]) = a(r)\V(r)\ 2 + ^(r)) 4 + 7 (r) \d^(v)\ 2 + -L(B(r) - H) 2 . (2.2) 

Here ^(r) and A(r) are the order parameter and 2D vector potential respectively, 

dA y dA x \ 



B = z 

and d- is the gauge invariant gradient 



dx dy 
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cL = -ih-^- + — A. (2.3) 
or c 

The space dependent Ginzburg-Landau coefficients a(r) and 7(1") for defects placed at the points 
{r,} have the form 
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7 (r) = 7o + ^ 7l (r-r i ), 
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where ao < and 70 > correspond to a uniform superconductor and the short-range functions 
ai(r) and 71 (r) describe the perturbation of these coefficients caused by a columnar defect located 
at the origin. 

According to the standard procedure one has to minimize the density (|2.2[) (i.e. to solve the 
Ginzburg-Landau equations), to find the extremal order parameter \1/ and vector potential A and 
to substitute them into the formula ( |2.1|) . Assume now that the density n of defects is small, 
n£ 2 (T) <C 1, i.e. the average distance between defects rC 1 ! 2 is much larger than the coherence 
length £(T) (which is of order of the distance a between neighboring vortices of the Abrikosov 
lattice). In this case the concentration expansionEl of the thermodynamic potential density in linear 
approximation yields 

F = f A + nJ{h{r)-f A )dv. (2.4) 

Here 

A _ 1 {H-H c2 f 
1 ~ 87r(2 K 2 -l)/V [ ^ 



is the free energy of the Abrikosov triangular lattice, Pa = 1.1598i3 and /i(r) is the minimum of the 
density (p.2|) containing a single defect placed at the origin (a(r) = ao + ai(r), 7 (r) = 70 + 7i(r)). 
Near the upper critical field 
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corresponding to a uniform superconductor with a(r) = ao, 7 (r) = 70, the minimization procedure 
can be applied to the density of the thermodynamic potential 

h (r; [a, 7, *]) = ao (1 - h) |*| 2 + | (l - ^5) |*| 4 + a x (r)|*| 2 + 7 i(r)|<9° *| 2 , (2.6) 
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which depends only on the order parameter!!!'!. Here h = H/H C 2, K is the Ginzburg-Landau param- 
eter 




and d°_ is defined by Eq.( [2.3| ) where the vector potential A of an applied field H stands for A. In 
what follows we will use the vector potential in the symmetric gauge A = H(—y/2, x/2). 

To find the order parameter which realizes this minimum one can use an expansion of \I/ in terms 
of Landau functions Lk(x, y) ( |A3| ) of the lowest Landau level of a particle with electron mass and the 
charge — 2e in the magnetic field H, substitute this expansion into Eq.( |2.6| ) and find the expansion 
coefficients from the minimum conditions. Such an expansion serves as a good approximation and 
one can neglect the contribution of the highest Landau levels even at a field H = 0.5# c2 §. In linear 
concentration approximation the problem is reduced to a single defect problem. Therefore in the case 
of isotropic functions «i(r) and 71 (r) the symmetry of the unperturbed Abrikosov lattice enables us 
to consider only two cases corresponding either to C§ symmetry, or to C3 symmetry. The hexagonal 
symmetry corresponds to the distorted vortex lattice with one vortex placed on the defect. The 
trigonal one corresponds to the lattice with the defect located in the center of the vortex triangle. 
In the hexagonal case the trial order parameter can be written as 

00 

* 6 (x, y)=i£ [vr- 1 M 6 (fc) + D{k)] L k (x, y). (2.7) 

k=0 

Here D(k) are the variational parameters which should be found. The case when all D are equal 
to zero and only the coefficients M§(k) remain, corresponds to the order parameter ty§(x,y) which 
describes the Abrikosov lattice with one of the vortices located at the origin and one of the symmetry 
axes parallel to the x-axis. The coefficients M§(k) (see Eq . (|A4j) of Appendix) are real and obey the 
selection rule! k = 6K +1, K = 0,1,2,... . In the trigonal case the trial order parameter is 
written as 

00 

y) = J2 r " [^Msik) + D(k)] L k (x, y). (2.8) 

k=0 

The case when all D are equal to zero, corresponds to the order parameter ^f(x, y) which describes 
the Abrikosov lattice whose origin coincides with the center of the vortex triangle and one of the 
symmetry axes is parallel to the x-axis. The real coefficients M 3 (k) ( |A5| ) obey the selection rule 
k = 3K,K = 0,1,2,... . 

Thus to obtain the lattice deformation caused by a single defect we have to find separately the 
extremal set of the variational parameters D(k) within each of the two symmetry classes separately, 
and to choose the most preferable one from the two of them. This procedure and its consequences 
will be discussed in the next two sections. 



III. WEAK DEFECTS 

Consider first a system with weak defects (in a sense that will be clear later on). It is natural 
to assume that in this case the two last terms in the thermodynamic potential density (|2.6|) do not 
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contribute to the variational equation for the order parameter and the latter one coincides with its 
Abrikosov value ^ A . Accordingly, the equilibrium thermodynamic potential can be written as 

F = f A + nJ ( ai (r)|^ A | 2 + 7i(r)|«9° V A \ 2 )dr. (3.1) 

Let us then specify the functions ai(r) and 71 (r) which describe the perturbation of the Ginzburg- 
Landau coefficients by defects 



ai{r) = — a aexp 



71 (r) = 7o7exp 
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Here a and 7 describe the strengths of the defect, measured in units ao and 70 respectively, and L 
is its size. Accurate estimations show that if the properly scaled strengths of defects 
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are small, a, 7 <C 1, then one can indeed neglect the Abrikosov lattice distortion. 

It seems that in the attractive case the preferable configuration is always a Cq, i.e. Abrikosov 
lattice with one of the vortices located on a defect. Nevertheless we will show that even in the case 
of small deformation, the previous statement is not always valid. In the general case one should take 
into account the two possible types of lattice symmetry with respect to a given defect, i.e. the two 
Abrikosov order parameters ty A 3 corresponding to the Cq symmetry 
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and to the C3 symmetry 



*i(x,y) =n- 1 ^2^ k M 3 (k)L k (x,y). 

k=0 



Substitution of these order parameters in the Eq. (|3.1|) yields the corresponding thermodynamic 
potentials, 

F* = f A + trn f A M 6 2 (ih^ 
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Fs = f A + ¥Ii f A Mi(0)a<p. (3.3) 
Here cp is the dimensionless scaled defect size 
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c = na 2 v3/2 is the dimensionless defect concentration (the number of defects per a single vortex) 
and a is the Abrikosov lattice constant (see Eq. ([Al| ) below). 

Suppose we deal with a defect in which the only variation parameter is a (i.e. the transition 
temperature). Then if a > such defect increases the thermodynamic potential F 3 leaving the _F 6 
unchanged ( |3.3| ). This means that the defect attracts a vortex and the symmetry Cq is preferable. 
Obviously, in the opposite case a < the symmetry C3 is preferable and therefore the defect 
is repulsive. This confirms the qualitative speculations presented in the Introduction. But the 
question is what happens if variation of both a and 7 is allowed. To answer this question we must 
compare the correction terms in Eqs.( [3.3|) . Taking values of Mq(1) and M^O) from the Tables I, II 
of Appendix we conclude that if 

7 < 0.387a, (3.5) 

then the symmetry Cq is preferable (attractive case). Yet, the ratio a/7 grows with the magnetic 
field (see Eq . (|3 . | ) above). Therefore if, for some field, 7 is slightly less than 0.387a then further 
increasing of a magnetic field can violate the inequality (|3.5|) and causes a first order phase transition 
to the C3 symmetry. 

In the region of field considered above, the vortex lattice is comparatively rigid and vortex repul- 
sion dominates above vortex-defect interaction. However even in this region the type of the lattice 
symmetry can be changed. For stronger magnetic fields or for stronger defects the lattice deforma- 
tion near defect is not negligible any more. This leads to richer and more complicated properties of 
the vortex system even in the case when 7 = 0. 



IV. STRONG DEFORMATION OF THE VORTEX LATTICE 



Now consider the case when a lattice deformation near defects is essential. This deformation is 
completely described by an infinite set of variational parameters {D(k)}. Direct substitution of the 
test function \l/(x, y) expressed in the forms (|2.7|) or Q2.8|) into the expression for the thermodynamic 
potential fl2j), (gT§) yields 



where 
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J(k + l)(D(k)D(l) + c.c) \ + ]T \D(k)\ 2 + 



^K'M^ + Dik)] 2 - 



+ <P) 1 



a<p + -/(ip 2 + k(l + 2<p 2 )) 
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and 



i(M) = E 



(/ + m)\ 



M(m)M(l + m- k) 



^Jk\l\m\{l + m-k)\ 



m = e 



k\ M(m)M(k - m) 
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This expression for the correction to the thermodynamic potential is general and valid for both two 
symmetries C 6 and C 3 . In each of these cases one should take into account the selection rules 



M 3 {k) = 5 k;3K M 3 {3K), 
M 6 (k) = 5 k , 6K+1 M 6 (6K + l) 
if = 0,1,2,... , 



(4.3) 



and use for M 3j6 (/c) their corresponding (real) values (see Eqs.([X4l), ( |A5| ) below). The next step is 
the minimization of the thermodynamic potential ( (4.1| ), ( |4.2| ) with respect to the coefficients {D(k)}. 
The equations which determine {D(k)} have the form 
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E (fc + m)!M(fc ~ / + m) D(l)D*(m) + E 
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[7iD(l)D(m)D*(l + m-k) + M(l + m- k)D{l)D{m)) 



2I(k,l)D(l) 



'(k + l)\ 
k\l\ 



J(k + l)D*(l) 



-D(k) + U- 1 M(k)+D{k) 



ao9 + 7(o9 2 + A;(l + 2a9 2 )) 







(4.4) 



and were obtained by Ovchinnikovi who used their linearized version for studying possible structural 
transitions. 



We numerically solve the infinite nonlinear system of Ovchinnikov equations without any simplifi- 
cation. The only (quite natural and verified) assumption which we use is that the perturbed lattice 
conserves its initial symmetry. This means that the coefficients {D(k)} obey the same selection rules 

D 6 (k) = 6 kt6K+1 D a (GK + l), 
D 3 {k) = 5 kj3K D 3 (3K), 



K 



0,1,2, 



that the initial coefficients M(k) do. Our strategy is as follows. For fixed values of the parameters 
a and 7 and for a fixed magnetic field we calculate the coefficients D(k) for two possible symmetries 
C 6 and C 3 . Then we substitute these solutions together with the corresponding sets of {M} into 
Eq.(O) and choose the most preferable solution which determine the vortex lattice deformation as 
well as the thermodynamics of the system to first order in the low concentration approximation. 
Thus to understand the results obtained we should first analyze the behavior of the coefficients D(k) 
in a magnetic field and to explain how this behavior influences to the order parameter evolution 
within each of the two symmetries separately. Then we can describe the vortex configuration, cor- 
responding to the preferable solution for a fixed set of parameters a, 7, L, and its evolution in a 
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magnetic field. 



The qualitative information concerning the behavior of the coefficients D(k) in a magnetic field can 
be obtained directly from Eqs.( |4.4j ). Consider for example an attractive defect with a > and 7 = 0. 
In this case, if one is not too close to the critical field H c2 the hexagonal symmetry should be realized 
and one starts from an analysis of the Cq solutions. Due to selection rules, the first nonvanishing 
equation of the system ( |4.4j ) will correspond to the value k = 1. This equation strongly depends 
on the (scaled) defect parameters a, 7 and ip, which are collected in the last term of Eq.fJOj). But 
right in the next equation (which corresponds to the value k = 7) this term is proportional to <p 7 
and due to the short range nature of the defect (ip < 1) is very small. Therefore all the higher order 
equations ( fOD with k = 13,19,... are practically homogeneous. As a result, the solution of fl4.4|) 
will give nonzero coefficients D(k) only for some small values of k. Thus the deformation of a vortex 
lattice happens mainly near the defect, at the distance of order of the Larmor radius Rk oc \fk~ Tl 



v max 



corresponding to the largest value of k such that D(k max ) 0, while the rest of the lattice remains 
undistorted. 
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FIG. 1. The dimensionless square modulus of the order parameter A in the hexagonal case for parameters 
a = 0.5, 7 = 0.0, tp/h = 0.5, h = 0.93. Seven vortices collapse on the defect. 

With raising of the applied magnetic field the effective coupling constants a and 7 increase dras- 
tically (see Eq. (|3.2| )), while the parameter <p fl3.4|) does not undergo any visible change. This leads 
to increasing values of the higher coefficients D(k) in the expansion ( |2.7|) of the order parameter 
and as a result, to spreading of the deformation far from the defect. The further the growth of the 
magnetic field is, the larger are the effective coupling constants. This implies that the last term in 
the Eq.( [4.4| ) for k = 1 becomes much larger than all preceding terms. In this case the solution is 
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Dq(1) = — 7r _1 M6(l), i.e. the first expansion coefficient practically reaches its limiting value. This 
value completely compensates the contribution of the unperturbed Abrikosov lattice to the k = 1 
expansion coefficient in Eq. ( |2.7|) . In this region of fields the expansion (|2.7|) begins from k = 7. The 
order parameter in the nearest vicinity of the defect becomes 



^ oc r 7 e 7vd . 

This means that the six nearest vortices have (almost) collapsed on the defect which pins the 
vortex containing seven flux quanta. One can see this effect on fig.|l| Here the quantity 

|A 6 | ~^cf, (4-5) 

which is proportional to the square modulus of the order parameter (normalization constant C is 
defined by Eq.(|A2|)), is plotted. 



With the further increasing the applied field the next coefficients D 6 (7), _D 6 (13), and so on will 
reach their limiting compensation values — 7r _1 M 6 (7), — 7r _1 M 6 (13), ... , and one could principally 
get a vortex containing thirteen, nineteen an so on flux quanta. However, numerical calculations 
show that for a realistic field range (not extremely close to the upper critical field) only the first 
collapse can be realized. 




x/a 

FIG. 2. The dimensionless square modulus of the order parameter A in the trigonal case for parameters 
a = 0.5, 7 = 0.0, ip/h = 0.5, h = 0.85. Three vortices collapse on the defect. 

A similar behavior of the expansion coefficients {D(k)} takes place in the trigonal case C3. Here 
in the case of attraction the coefficient D 3 (0) is the first one which reaches its compensation value 







— n^ 1 M^(0), that corresponds to the three vortices collapse on the defect. Such a configuration is 
displayed on fig.^ where the quantity |A 3 | 2 defined by the r.h.s of Eq.( [4.5|) , with ^/q replaced by ^3, is 
plotted for the same values of parameters as in the hexagonal case and for the applied field h « 0.85. 
With increasing of the magnetic field one expects the appearance of six-, and so on multy-quanta 
vortices. As in the previous case, numerical analysis shows that only the first collapse happens in a 
realistic range of field. 
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FIG. 3. The dimensionless square modulus of the order parameter A in the hexagonal case for parameters 
a = 0.5, 7 = 0.0, (f/h = 0.5, h = 0.85. Attractive defect causes a comparatively weak vortex lattice deformation. 

Note that for the same set of parameters the first collapse within the trigonal symmetry occurs at 
a weaker field (h « 0.85) than in the hexagonal symmetry (h m 0.93). The reason is that in the Cq 
system seven vortices must overcome their mutual repulsion in order to fall on the defect, while in 
the C3 system only three vortices collapse. For the field h ~ 0.85, at which, in the symmetry C3, 
three vortices are already collapsed on the defect (fig.@), in the Cq symmetry, the lattice is distorted 
but still without any vortex collapse (fig||). 

The vortex lattice deformation near a repulsive defect is presented in fig.|j. Here the three nearest 
vortices are slightly shifted from the defect and a visible deformation occurs only in the nearest 
vicinity of the defect. 

Up to now we analyzed the solutions of Eqs.( |4.4| ) within two symmetries Cq and C3 separately. 
Now we can choose the most preferable one from them and describe the typical vortex lattice behav- 
ior in some interval of the magnetic fields close to the upper critical field. We start from the same 
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case of attractive defects a > (7=0) of a small concentration. If an applied field is not too close to 
H c2 , then a deformation of the lattice near a single deffect is small and the preferable local symmetry 
near each defect is C 6 . The defects are occupied by vortices and the rest of the lattice is slightly 
deformed. With increasing of the magnetic field the deformation near defects becomes stronger (as 
shown in fig.|3]) and at some critical field hi the C3 solution of Eqs. ( |4.4| ) corresponding to collapse of 
three vortices on the defect becomes preferable (see fig.0)- As a result, a local structural transition 
Cq — > C3 occurs. With further increasing of the field, one deals with C3 symmetry, three vortices 
occupying the defect and the deformation of the nearest (with respect to the defect) part of the 
vortex lattice is observed. But at some critical field /12 the Cq solution of Eqs. ( |4.4| ) corresponding to 
collapse of seven vortices on the defect (see fig.[I|) becomes preferable and a local structural transi- 
tion C3 — > Cq occurs and so on. Thus, one has a sequence of reentering first order phase transitions 
Cq ~~ > C3 ~~ > C*6 ~~ > ■■■ ■ 
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FIG. 4. The dimensionless square modulus of the order parameter A in the hexagonal repulsive case for parameters 
a = -0.5,7 = 0.0, (p = 0.5. 

A similar analysis can be done in the general case where q^O and 7 7^ 0. The numerical results 
obtained for various sets of parameters and magnetic field enables us to construct a phase diagram in 
the (a, 7) plane for a fixed scaled size ip (|3.4j ) of a defect. Part of such diagram is given in fig.|5|. Here 
the two solid curves separate the regions where the local symmetry is hexagonal (Cq) or trigonal 
(C3) . Near the upper critical field ip ~ L 2 /£ 2 and the diagram becomes universal. For each fixed 
defect parameters and for each value of the magnetic field the diagram enable us to determine the 
preferable local symmetry of the system. 

To explain how to extract this information from the phase diagram consider a sample with some 
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fixed parameters a, 7, and L, and start from an initial applied field h = H /H C 2. This corresponds 
to a starting point (a , 70) in the diagram of fig.|| where a an d 7o are determined by Eqs. (|3.2| ) with 
h = h . Further evolution of the parameters a and 7 with growth of the magnetic field is described 
by equation 

7 

7 = - (" - «o) + 7o 

a 

and corresponds to some ray on the phase diagram, starting at the initial point (a , 70) and directed 
out of the origin. Four such rays are displayed in fig.|5]. For all rays the starting field is ho = 0.9 and 
a = 0.1. The increasing of the magnetic field leads to the change in the effective coupling constants 
( p.2|) i.e. to the motion of a starting point along the ray. This movement in its turn results in a 
sequence of reentering transitions from one local symmetry to another. 




FIG. 5. Phase diagram of superconductor on (0,7) plane. All rays: a — 0.1; ho = 0.9. Solid ray 7 = 0.01. Dashed ray 
7 = 0.03. Dotted ray 7 = 0.06. Dashed-dotted ray 7 = -0.01. 

The solid ray corresponds to 7 = 0.01. It is seen from fig.^ that at the initial field ho = 0.9 the local 
deformation has a C§ symmetry. This is consistent with the analytical prediction of Section III: the 
deformation around a defect is small and the inequality ( p.5| ) is valid. As the field increases, the ray 
crosses the lower solid curve and the sample undergoes a first order phase transition to the trigonal 
local symmetry C3. At this symmetry we have three collapsed vortices at each defect. Transition 
to the symmetry C$ back is also possible, but it is not seen on the diagram because it occurs in the 
region a>l i.e. at a field extremely close to the upper critical field H c2 . 

The dashed ray corresponds to the value 7 = 0.03 and represents probably the most interesting 
case. Here even in the comparatively low field h « 0.775 < ho (the corresponding point of the ray is 
not displayed on fig.^D the Cq-C^ symmetry transition occurs. In both the two lattice configurations 
below and above the transition the lattice deformation is small and can be described within the 
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approach of Section III. The inequality ( |3.5|) is violated below the transition field and is valid above 
it. The dashed ray on the diagram starts from the field h Q = 0.9 and for the first time crosses the 
lower solid curve at a field h « 0.906, at which the lattice undergoes the next C3 — > C$ transition. No 
vortex collapse still happens at this field because the value of Dq(1) is still far from its compensating 
value. However two next transitions take place because of vortex collapse. The second transition to 
the symmetry C3 at a field h m 0.94 happens when the coefficient D^(0) in the symmetry C3 almost 
reaches its compensating value 1)3(0) = —ir M^O) and therefore this transition corresponds to the 
collapse of the three vortices at the defect. Similarly the third transition to the symmetry Cq at a 
field h ~ 0.99 corresponds to the collapse of the seven vortices at the defect. 

The dotted ray in fig.||] corresponds to 7 = 0.06 that provides only a C3 symmetry in the compar- 
atively low field region. In the high field region we obtain C3 — > Cg transition at the field h ~ 0.98. 

Note that the figures || and || already reffered to above, present the contour plots of the order 
parameter near defect in the vicinity of the Cq — > C3 transition due to collapse of the three nearest 
vortices on the defect. These plots correspond to the point (« 4, 0) on the ray coinciding with 
the positive a— semiaxis on the phase diagram. At this point the order parameter exhibits a small 
deformation in the symmetry Cq as it is displayed in fig.0, while in the symmetry C3 it is strongly 
deformed due to the collapse (see fig. Q). 

In the region where a > and 7 < a local symmetry transition Cg — > C3 due to vortex collapse 
is described by the last fourth ray on the diagram. This ray corresponds to the parameters a = 0.1, 
7 = -0.01. 



V. SMALL CONCENTRATION OF THE DEFECTS 

During the two previous sections we dealt with a single defect problem. To be sure that our re- 
sults Q4.ip do describe a macroscopic system with a finite concentration of the defects we have to be 
sure that the next (second order) concentration correction to the thermodynamic potential is small. 
To estimate this correction one has to solve exactly the two defects problem which is much more 
complicated. Therefore we choose another way. 

Consider for simplicity an attractive case and magnetic field which is not too close to H c2 . Put 
the undistorted vortex lattice on the plane where (point) defects are distributed and shift one of 
the vortices nearest to each inhomogeneity to the position of that inhomogeneity. There are many 
similar ways to arrange the vortex lattice, but one has to choose such a way which leads to alter- 
nation of the regions where the lattice is compressed with ones where its rarefied. Finally let us 
distort the regions of the lattice close to inhomogeneities according to the results obtained within 
single defect approximation. This latter distortion is already taken into account exactly. So one has 
only to estimate the additional contribution to the thermodynamic potential from the intermediate 
regions (between inhomogeneities) whose deformation is well described by elastic theory. 

The number of extra vortices per region is of order of unity. Therefore the deformation tensorup to 
a numerical factor of the order of unity equals to dimensionless concentration of defects c. The cor- 
rection to the thermodynamic potential will be of the order of Cc 2 , where C is the elastic modulus. 
But the elastic part of the deformation has an alternating behavior with a characteristic wavelength 
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of the order of the average distance between inhomogeneities. As it was shown by E. Brandtfl, all 
the elastic moduli are proportional to (1 — H/H c2 ) 2 if this distance is much less than the penetration 
length divided by (1 — H/ H c2 Y^ 2 . The latter inequality can be rewritten as 36kc ^> 1 — H/H c2 , 
where k 3> 1 is the Ginzburg-Landau parameter. In the region of parameters which we are mostly 
interested in c = 0.03, 1 — H/H c2 = 0.06 and the inequality k ^> 1 is evidently valid. This means 
that corresponding contribution to the thermodynamic potential is of the order of (1 — H/H C 2) 2 c 2 . 
This is exactly the second order concentration correction which in the case c 1 is smaller than 
the contribution accounted for within the linear concentration expansion. 

Thus in the case of small concentration one can use the results obtained in the two previous sections 
and describe the thermodynamics of the system near H c2 . Define a dimensionless magnetization 

m = -Au{2k 2 - 1){3a^~ 

and dimensionless magnetic succeptibility 

dm 




All the local symmetry transitions described above manifest themselves as jumps on the magne- 
tization curve (figJI) and as discontinuities on the magnetic succeptibility curve (fig.0). The most 
pronounced jumps occur at the two transitions accompanied by vortex collapse, namely at the fields 
h w 0.94 and h w 0.99. 



14 



O.O ts» 



-0.2- 



-0.4- 




X 



-0.6 - 



-0.8 - 



-l.O- 



0.90 0.92 0.94 0.96 0.98 1 .OO 

H/H 



FIG. 7. The magnetic succeptibility of superconductor for the parameters alpha = 0.1,7 = 0.03 and for the concentration 
values c = 0.03 (triangles) and c = 0.05 (circles). 



VI. SUMMARY 



We studied the equilibrium properties of conventional type II superconductor with small concen- 
tration of randomly placed identical columnar defects. In the vicinity of the upper critical field the 
vortex lattice undergoes a strong deformation with two possible local symmetries - hexagonal one 
Cq and trigonal one C3. The character of the deformation is determined by the vortex-defect inter- 
action. The vortices can collapse onto attractive defects and the formation of multy-quanta vortices 
becomes possible. Formation of the multiquantum vortices was predicted earlier^, but in "twice" 
opposite limiting case. We deal with a short-range defect and gain an energy because of softening 
of the Abrikosov lattice near H c2 , while int3 a very strong defect with a radius comparable with the 
penetration length was considered. 



Increasing the external field gives rise to the reentering transitions between the two possible types 
of symmetry. These transitions can be described by a universal phase diagram. They manifest 
themselves as jumps of the magnetization and peculiarities of the magnetic susceptibility 

One of the way to observe these equilibrium states near H c2 is to cool a sample subject to a 
magnetic field in the normal state, below the critical temperature. Another possibility is to observe 
not the equilibrium state as a whole, but visualize the local deformation of the vortex lattice near 
defects. 

VII. ACKNOWLEDGEMENTS 



This study was partially supported by the Israel Science Foundation. We are grateful to H. Brandt, 
E. Chudnovsky, B. Horovitz, V. Mineev, R. Mints, Z. Ovadyahu, B.Ya. Shapiro, and E. Zeldov for 
very helpful discussions of our results. 



15 



APPENDIX A: ABRIKOSOV LATTICE EXPANSION 



The Abrikosov order parametera 

^I A (x,y) = Cexp 



03 I % - 

a a A 



which was obtained in the Landau gauge A = xH(—y, 0) describes a triangular vortex lattice with 
sites 

j x = a(2m + n + 2)/2, 
1 y = aV3(2n + l)/4, 



(m and n are integers). Here 



a(T) = 2Z(T). 



7T 



h(T)y/3 



(Al) 



is the triangle side (Abrikosov lattice constant), and /i(T) = H/H C 2(T). The Abrikosov normalization 
constant C is related to the thermodynamic potential density ( |2.5|) of the clean superconductor by 



|C| 



2 3^21/^1 



|ao|(l - fc) 



(A2) 



and 9z(u\t) is the Euler 9 -functionEl 
The order parameter 



.i'o.;yo 



>, y) = exp 



2tt .4tt?/o 
— t= — xy + i — 7^X 



y/3t 



a* 



y/3t 



a 



^/ A (x - ax ,y- ay Q ), 



describes the shifted Abrikosov lattice in the symmetric gauge. This order parameter can be ex- 
panded as 



k=0 



with respect to Landau functions with the orbital moment k > 



Lk(x,y) 



C 2vr 



3V4 V jfe! 



\/2~k r 


k 


. 31/4 a . 


exp 



(A3) 



of the lowest Landau level of a particle with electron mass and charge — 2e in the magnetic field H. 
The expansion coefficients are 



M xo ,y (.k) 



k~2m 



k—2m 



n=— 00 m=0 



.7T 



2 m m\(k - 2m)! 



n + 



k—2m 



X 



n — ri(l + 4z. 



0, 



TTV 7 ^ / 2?/ N 
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The case x = —1/2, y = —a/3/4 corresponds to the C 6 symmetry when one vortex is placed at 
the origin. Taking into account the selection rule ( [4.3| ) we write down the expansion coefficients as 



where 



M 6 (k) 



oo fc-2m>0 % 

E E~ 



fc-2m-l 



2ny/3 



k—2m 



n=—oo m=0 



2 m m\(k - 2m)! 



n 



k—2m 



exp i\ 



n + n 



71 



n 

2 



(A4) 



The expansion coefficients M 6 (/c) are real. The values of the first few of them are contained in the 
table I. 



TABLE I. Values of an Expansion Coefficients Mg(fc) 



k M 6 (fc) 

1 -2.792 

7 4.057 

13 2.260 

19 -4.852 

25 -1.494 

31 -3.538 

37 4.817 

43 2.605 

49 0.479 

55 5.180 



The second case xq = —1/2 and yo = \/3/12 corresponds to the C3 symmetry, when the origin of 
the coordinate system is placed in the center of an elementary triangle. Selection rules allow us to 
write down the expansion coefficients as 

M_ 1/2 ^- 3/12 (k) = t- 3K 5 k , 3K M 3 (3K). 
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The explicit expression for M 3 (k) is given by the formula 



M 3 (k) 



fc-2m>0 I 



2fc-2m 



E E 

n=— oo m=0 



27^ 



k—2m 



2 m m\(k - 2m)! 



n + 



6 



k—2m 



X 



7T 



exp < 



n + n 



(AS) 



All the coefficients M 3 (k) also are real. The values of the first few coefficients are given by Table II. 

TABLE II. Values of an Expansion Coefficients Ms(k) 



k_ 


3 

6 

9 

12 

15 

18 

21 

24 

27 

30 



M 3 (k) 

1.738 

-2.942 

-2.227 

1.646 

-3.568 

-2.310 

-0.756 

3.185 

-1.563 

3.110 

3.606 



1 L.V. Shubnikov, V.I. Khotkevich, Yu.D. Shepelev, Yu.N. Ryabinin, Zh. Eksper. Teor. Fiz. 7, 221 (1937). 

2 A. A. Abrikosov, Sov. Phys. JETP 5, 1174 (1957). 

3 A.I. Larkin, Yu.N. Ovchinnikov, in Nonequilibrium Superconductivity, Eds. D.N. Landerberg and A.I. Larkin (Elsevier, 
Amsterdam, 1986). 

4 A.I. Larkin and Yu.N. Ovchinnikov, J. Low Temp. Phys. 34, 409 (1978). 

5 R. Labusch, Phys. Stat. Sol. 32, 439 (1969). 

6 A.I. Larkin, Sov. Phys. JETP 31, 784 (1970). 

7 E. H. Brandt, J. Low Temp. Phys. 26, 709, 735 (1977). 

8 Yu.N. Ovchinnikov, Sov. Phys. JETP 55, 1162 (1982). 

9 J.G. Bednortz, K.A. Muller, Z. Phys. 64, 189 (1986). 

G. Blatter,M.V. Feigelman, V.B. Geshkenbein, A.I. Larkin, V.M. Vinokur, Rev. Mod. Phys. 66,1125 (1994). 

1 E.H. Brandt, Rep. Progr. Phys. 58, 1465 (1995) 

2 a) L. Civale, A.D.Marwick, M.W. McElfresh, T.K. Worthington, A. P. Malozemoff, F.H. Holtzberg, J.R. Thompson, M.A. 
Kirk, Phys. Rev. Lett. 65, 1164 (1990); b) L. Civale, A.D.Marwick, T.K. Worthington, M.A. Kirk, J.R. Thompson, L. 
Krusin-Elbaum, Y. Sun, J.R. Clem, F.H. Holtzberg, Phys. Rev. Lett. 67, 648 (1991). 

3 A.I. Buzdin, Phys. Rev. B 47 11416 (1993) 

4 M. Baert, V.V. Metlushko, R. Jonckheere, V.V. Moshchalkov, Y. Bruynseraede, Phys. Rev. Lett. 74, 3269 (1995). 

5 LB. Khalfin, BYa. Shapiro, Physica C 202, 393 (1992). 

6 V.L. Ginzburg, L.D. Landau, Zh. Eksper. Teor. Fiz. 20, 1064 (1950). 
7 I.M. Lifshits, Nuovo Cimento 3, Suppl, 716 (1963). 

8 W.M. Kleiner, L.M. Roth, S.H. Autler, Phys. Rev. A133, 1226 (1964). 

9 E.H. Brandt, Phys. Stat. Sol. (b) 71, 277 (1975). 

10 E.H. Brandt, Phys. Stat. Sol. (b) 51, 345 (1972). 

11 I.S. Gradstein, and I.M. Ryzhik, Tables of Integtrals, Sums, Series and Products (Academic, New York, 1980). 



18 



